Dərs 16: Graded Response Model (GRM)
Polytomous Item Response Theory Modellərinin Əsasları
IRT Kursu
2025-07-19
1 Dərsin Məqsədləri
Bu dərsin sonunda tələbələr:
✓ GRM mathematical formulation və polytomous IRT
prinsiplərini anlaya bilə
✓ Category Response Functions (CRFs) yarada və
interpretasiya edəcək
✓ Threshold parameters və onların mənaını izah
edəcək
✓ Rating Scale Model və Partial Credit
Model ilə GRM-i müqayisə edəcək
✓ R paketlərində (ltm, mirt, eRm) polytomous IRT
analysis aparacaq
✓ Model fit və model selection
metodlarını tətbiq edəcək
✓ Practical applications ilə real data-da GRM istifadə
edəcək
2 GRM-nin Riyazi İfadəsi
2.1 Graded Response Model Mathematical Foundation
## === GRM MATHEMATICAL FOUNDATION ===## Graded Response Model: Ordered categorical responses üçün IRT model## • Samejima (1969) tərəfindən yaradılmış## • Cumulative probability approach istifadə edir## • Discrimination və threshold parameterlər## • Likert scales və rating data üçün ideal# Comprehensive GRM mathematical demonstration
demonstrate_grm_mathematics <- function() {
cat("=== GRM MATHEMATICAL DEMONSTRASIYASI ===\n")
# Generate example GRM data
grm_data <- generate_grm_data(
n_items = 6,
n_persons = 1500,
n_categories = 5, # 0, 1, 2, 3, 4 (5-point Likert)
seed = 1601
)
cat("GRM Example Dataset:\n")
cat("Items:", ncol(grm_data$responses), "\n")
cat("Persons:", nrow(grm_data$responses), "\n")
cat("Categories per item:", grm_data$n_categories, "(0-4)\n")
cat("Response range:", range(grm_data$responses), "\n\n")
# Display mathematical formulation
cat("=== GRM MATHEMATICAL FORMULATION ===\n")
mathematical_formulation <- c(
"**CUMULATIVE PROBABILITY MODEL:**",
"",
"P*ik(θ) = P(Xij ≥ k | θi) = exp[ai(θi - bik)] / [1 + exp[ai(θi - bik)]]",
"",
"Where:",
"• P*ik(θ) = Cumulative probability of responding in category k or higher",
"• ai = Discrimination parameter for item i",
"• bik = Threshold parameter k for item i",
"• θi = Ability level of person i",
"",
"**CATEGORY PROBABILITIES:**",
"",
"P(Xij = k | θi) = P*ik(θ) - P*i(k+1)(θ)",
"",
"Special cases:",
"• P(Xij = 0 | θi) = 1 - P*i1(θ)",
"• P(Xij = mi | θi) = P*imi(θ) [highest category]",
"",
"**CONSTRAINTS:**",
"",
"• bi1 < bi2 < ... < bi(mi) [Ordered thresholds]",
"• Σk P(Xij = k | θi) = 1 [Probabilities sum to 1]",
"• All P(Xij = k | θi) ≥ 0 [Non-negative probabilities]"
)
for(line in mathematical_formulation) {
cat(line, "\n")
}
# Demonstrate with specific example
cat("\n=== NUMERICAL EXAMPLE ===\n")
# Use first item from generated data
example_item <- 1
a_example <- grm_data$discriminations[example_item]
b_example <- grm_data$thresholds[example_item, ]
cat("Example Item", example_item, ":\n")
cat("Discrimination (a):", round(a_example, 3), "\n")
cat("Thresholds (b):", paste(round(b_example, 3), collapse = ", "), "\n")
# Calculate probabilities at different theta levels
theta_examples <- c(-2, -1, 0, 1, 2)
probability_table <- data.frame()
for(theta in theta_examples) {
# Calculate cumulative probabilities
cum_probs <- numeric(grm_data$n_categories)
for(k in 1:(grm_data$n_categories - 1)) {
cum_probs[k] <- exp(a_example * (theta - b_example[k])) /
(1 + exp(a_example * (theta - b_example[k])))
}
cum_probs[grm_data$n_categories] <- 0
# Calculate category probabilities
cat_probs <- numeric(grm_data$n_categories)
cat_probs[1] <- 1 - cum_probs[1]
for(k in 2:(grm_data$n_categories - 1)) {
cat_probs[k] <- cum_probs[k-1] - cum_probs[k]
}
cat_probs[grm_data$n_categories] <- cum_probs[grm_data$n_categories - 1]
# Add to table
row_data <- data.frame(
Theta = theta,
P_Cat0 = round(cat_probs[1], 3),
P_Cat1 = round(cat_probs[2], 3),
P_Cat2 = round(cat_probs[3], 3),
P_Cat3 = round(cat_probs[4], 3),
P_Cat4 = round(cat_probs[5], 3)
)
probability_table <- rbind(probability_table, row_data)
}
cat("\nCategory Probabilities by Ability Level:\n")
print(probability_table)
# Verify probabilities sum to 1
row_sums <- rowSums(probability_table[, -1])
cat("\nRow sums (should all be 1.000):", paste(round(row_sums, 3), collapse = ", "), "\n")
# Expected score calculation
cat("\n=== EXPECTED SCORE FUNCTION ===\n")
expected_scores <- numeric(length(theta_examples))
for(i in 1:length(theta_examples)) {
theta <- theta_examples[i]
cat_probs <- grm_probability(theta, a_example, b_example)
expected_scores[i] <- sum(cat_probs * (0:(grm_data$n_categories - 1)))
}
expected_score_table <- data.frame(
Theta = theta_examples,
Expected_Score = round(expected_scores, 3),
Min_Possible = 0,
Max_Possible = grm_data$n_categories - 1
)
cat("Expected Score Function:\n")
print(expected_score_table)
cat("\nExpected Score Properties:\n")
cat("• Monotonically increasing with θ\n")
cat("• Bounded between 0 and", grm_data$n_categories - 1, "\n")
cat("• Smooth, continuous function\n")
cat("• Steepness determined by discrimination parameter\n")
return(list(
grm_data = grm_data,
example_item = example_item,
example_params = list(discrimination = a_example, thresholds = b_example),
probability_table = probability_table,
expected_score_table = expected_score_table
))
}
# Run GRM mathematical demonstration
grm_math_results <- demonstrate_grm_mathematics()## === GRM MATHEMATICAL DEMONSTRASIYASI ===
## GRM Example Dataset:
## Items: 6
## Persons: 1500
## Categories per item: 5 (0-4)
## Response range: 0 4
##
## === GRM MATHEMATICAL FORMULATION ===
## **CUMULATIVE PROBABILITY MODEL:**
##
## P*ik(θ) = P(Xij ≥ k | θi) = exp[ai(θi - bik)] / [1 + exp[ai(θi - bik)]]
##
## Where:
## • P*ik(θ) = Cumulative probability of responding in category k or higher
## • ai = Discrimination parameter for item i
## • bik = Threshold parameter k for item i
## • θi = Ability level of person i
##
## **CATEGORY PROBABILITIES:**
##
## P(Xij = k | θi) = P*ik(θ) - P*i(k+1)(θ)
##
## Special cases:
## • P(Xij = 0 | θi) = 1 - P*i1(θ)
## • P(Xij = mi | θi) = P*imi(θ) [highest category]
##
## **CONSTRAINTS:**
##
## • bi1 < bi2 < ... < bi(mi) [Ordered thresholds]
## • Σk P(Xij = k | θi) = 1 [Probabilities sum to 1]
## • All P(Xij = k | θi) ≥ 0 [Non-negative probabilities]
##
## === NUMERICAL EXAMPLE ===
## Example Item 1 :
## Discrimination (a): 0.811
## Thresholds (b): -0.935, -0.217, 0.068, 0.742
##
## Category Probabilities by Ability Level:
## Theta P_Cat0 P_Cat1 P_Cat2 P_Cat3 P_Cat4
## 1 -2 0.704 0.106 0.033 0.060 0.098
## 2 -1 0.513 0.141 0.050 0.100 0.196
## 3 0 0.319 0.137 0.057 0.132 0.354
## 4 1 0.172 0.099 0.048 0.128 0.552
## 5 2 0.085 0.057 0.030 0.092 0.735
##
## Row sums (should all be 1.000): 1.001, 1, 0.999, 0.999, 0.999
##
## === EXPECTED SCORE FUNCTION ===
## Expected Score Function:
## Theta Expected_Score Min_Possible Max_Possible
## 1 -2 0.742 0 4
## 2 -1 1.325 0 4
## 3 0 2.065 0 4
## 4 1 2.789 0 4
## 5 2 3.336 0 4
##
## Expected Score Properties:
## • Monotonically increasing with θ
## • Bounded between 0 and 4
## • Smooth, continuous function
## • Steepness determined by discrimination parameter##
## === GRM MODEL IDENTIFICATION ===# Model identification discussion
identification_constraints <- c(
"**IDENTIFICATION CONSTRAINTS IN GRM:**",
"",
"1. **Ability Scale:**",
" • Mean(θ) = 0, Var(θ) = 1 [Standard normal distribution]",
" • Or fix two anchor items",
"",
"2. **Threshold Ordering:**",
" • bi1 < bi2 < ... < bimi for each item i",
" • Ensures proper category ordering",
"",
"3. **Parameter Constraints:**",
" • Discrimination parameters ai > 0",
" • No upper bound on discrimination",
" • Threshold parameters unbounded but ordered",
"",
"**COMMON PARAMETERIZATIONS:**",
"",
"1. **Slope-Intercept (SI):**",
" • P*ik = exp(ai*θ + cik) / [1 + exp(ai*θ + cik)]",
" • Where cik = -ai*bik",
"",
"2. **Slope-Threshold (ST):**",
" • P*ik = exp[ai(θ - bik)] / [1 + exp[ai(θ - bik)]]",
" • Most common in IRT software",
"",
"**REPARAMETERIZATION RELATIONSHIPS:**",
"• SI to ST: bik = -cik/ai",
"• ST to SI: cik = -ai*bik"
)
for(constraint in identification_constraints) {
cat(constraint, "\n")
}## **IDENTIFICATION CONSTRAINTS IN GRM:**
##
## 1. **Ability Scale:**
## • Mean(θ) = 0, Var(θ) = 1 [Standard normal distribution]
## • Or fix two anchor items
##
## 2. **Threshold Ordering:**
## • bi1 < bi2 < ... < bimi for each item i
## • Ensures proper category ordering
##
## 3. **Parameter Constraints:**
## • Discrimination parameters ai > 0
## • No upper bound on discrimination
## • Threshold parameters unbounded but ordered
##
## **COMMON PARAMETERIZATIONS:**
##
## 1. **Slope-Intercept (SI):**
## • P*ik = exp(ai*θ + cik) / [1 + exp(ai*θ + cik)]
## • Where cik = -ai*bik
##
## 2. **Slope-Threshold (ST):**
## • P*ik = exp[ai(θ - bik)] / [1 + exp[ai(θ - bik)]]
## • Most common in IRT software
##
## **REPARAMETERIZATION RELATIONSHIPS:**
## • SI to ST: bik = -cik/ai
## • ST to SI: cik = -ai*bik##
## === GRM PARAMETER INTERPRETATION ===# Parameter interpretation guide
parameter_interpretation <- data.frame(
Parameter = c("ai (Discrimination)", "bik (Threshold k)", "Expected Score", "Information"),
Interpretation = c("Item's ability to differentiate between ability levels",
"Ability level where P(X≥k) = 0.5",
"Average response for given ability",
"Precision of ability estimation"),
Range = c("0 < ai < ∞ (typically 0.5-3.0)",
"-∞ < bik < ∞ (typically -3 to +3)",
"0 ≤ E[X|θ] ≤ m (max category)",
"0 ≤ I(θ) < ∞"),
High_Values = c("Steep expected score function",
"High ability needed for category k+",
"High ability level",
"Precise ability estimation"),
Low_Values = c("Flat expected score function",
"Low ability sufficient for category k+",
"Low ability level",
"Imprecise ability estimation")
)
print(parameter_interpretation)## Parameter Interpretation
## 1 ai (Discrimination) Item's ability to differentiate between ability levels
## 2 bik (Threshold k) Ability level where P(X≥k) = 0.5
## 3 Expected Score Average response for given ability
## 4 Information Precision of ability estimation
## Range High_Values
## 1 0 < ai < ∞ (typically 0.5-3.0) Steep expected score function
## 2 -∞ < bik < ∞ (typically -3 to +3) High ability needed for category k+
## 3 0 ≤ E[X|θ] ≤ m (max category) High ability level
## 4 0 ≤ I(θ) < ∞ Precise ability estimation
## Low_Values
## 1 Flat expected score function
## 2 Low ability sufficient for category k+
## 3 Low ability level
## 4 Imprecise ability estimation3 Category Response Functions
3.1 Understanding CRFs və Visualization
## === CATEGORY RESPONSE FUNCTIONS ===## Category Response Functions (CRFs): Hər kateqoriya üçün ehtimal əyriləri## • Müxtəlif ability səviyyələrində kateqoriya seçim ehtimalları## • Item charakteristikasının qrafik təmsili## • Threshold parameters və discrimination təsiri## • Model fit və item functioning qiymətləndirməsi# Comprehensive CRF demonstration
demonstrate_category_response_functions <- function() {
cat("=== CATEGORY RESPONSE FUNCTIONS DEMONSTRASIYASI ===\n")
# Create examples with different parameter configurations
theta_range <- seq(-4, 4, 0.1)
# Example 1: Well-functioning item
example1_params <- list(
discrimination = 1.5,
thresholds = c(-1.5, -0.5, 0.5, 1.5),
title = "Well-Functioning Item"
)
# Example 2: Low discrimination item
example2_params <- list(
discrimination = 0.8,
thresholds = c(-1.5, -0.5, 0.5, 1.5),
title = "Low Discrimination Item"
)
# Example 3: Highly discriminating item
example3_params <- list(
discrimination = 2.5,
thresholds = c(-1.5, -0.5, 0.5, 1.5),
title = "High Discrimination Item"
)
# Example 4: Poorly spaced thresholds
example4_params <- list(
discrimination = 1.5,
thresholds = c(-2.0, -1.8, -1.6, -1.4),
title = "Poorly Spaced Thresholds"
)
examples <- list(example1_params, example2_params, example3_params, example4_params)
# Calculate CRFs for each example
calculate_crfs <- function(params, theta_range) {
crf_data <- data.frame()
for(theta in theta_range) {
cat_probs <- grm_probability(theta, params$discrimination, params$thresholds)
for(cat in 1:length(cat_probs)) {
crf_data <- rbind(crf_data, data.frame(
Theta = theta,
Category = paste0("Cat", cat - 1),
Probability = cat_probs[cat],
Item = params$title
))
}
}
return(crf_data)
}
# Generate CRF data for all examples
all_crf_data <- data.frame()
for(i in 1:length(examples)) {
example_crf <- calculate_crfs(examples[[i]], theta_range)
all_crf_data <- rbind(all_crf_data, example_crf)
}
# Create CRF plots
if(require(ggplot2, quietly = TRUE)) {
cat("Creating Category Response Function plots...\n")
# Plot for each example
for(i in 1:length(examples)) {
example_data <- all_crf_data[all_crf_data$Item == examples[[i]]$title, ]
p <- ggplot(example_data, aes(x = Theta, y = Probability, color = Category)) +
geom_line(size = 1.2) +
labs(title = paste("Category Response Functions:", examples[[i]]$title),
subtitle = paste("Discrimination =", examples[[i]]$discrimination,
"| Thresholds =", paste(examples[[i]]$thresholds, collapse = ", ")),
x = "Ability (θ)",
y = "Probability",
color = "Category") +
scale_color_brewer(type = "qual", palette = "Set1") +
theme_minimal() +
theme(plot.title = element_text(hjust = 0.5, size = 14, face = "bold"),
plot.subtitle = element_text(hjust = 0.5, size = 11),
legend.position = "bottom") +
ylim(0, 1) +
xlim(-4, 4)
print(p)
}
}
# Analyze CRF characteristics
cat("\n=== CRF CHARACTERISTICS ANALYSIS ===\n")
analyze_crf_properties <- function(params, theta_range) {
# Calculate modal categories across theta range
modal_categories <- numeric(length(theta_range))
max_probabilities <- numeric(length(theta_range))
for(i in 1:length(theta_range)) {
cat_probs <- grm_probability(theta_range[i], params$discrimination, params$thresholds)
modal_categories[i] <- which.max(cat_probs) - 1 # 0-indexed
max_probabilities[i] <- max(cat_probs)
}
# Find threshold crossing points (where adjacent categories have equal probability)
crossing_points <- numeric(length(params$thresholds))
for(k in 1:length(params$thresholds)) {
# Find theta where P(X=k-1) = P(X=k)
theta_k <- params$thresholds[k] # Approximate crossing point
# Refine with optimization
objective <- function(theta) {
cat_probs <- grm_probability(theta, params$discrimination, params$thresholds)
abs(cat_probs[k] - cat_probs[k+1])
}
tryCatch({
result <- optimize(objective, interval = c(theta_k - 1, theta_k + 1))
crossing_points[k] <- result$minimum
}, error = function(e) {
crossing_points[k] <- theta_k
})
}
# Category usage statistics
category_usage <- table(modal_categories)
n_categories_used <- length(category_usage)
properties <- list(
modal_categories = modal_categories,
max_probabilities = max_probabilities,
crossing_points = crossing_points,
category_usage = category_usage,
n_categories_used = n_categories_used,
min_max_probability = min(max_probabilities),
effective_range = range(theta_range[max_probabilities > 0.4])
)
return(properties)
}
# Analyze each example
for(i in 1:length(examples)) {
cat("Analysis for", examples[[i]]$title, ":\n")
properties <- analyze_crf_properties(examples[[i]], theta_range)
cat("Categories used:", paste(names(properties$category_usage), collapse = ", "), "\n")
cat("Crossing points:", paste(round(properties$crossing_points, 2), collapse = ", "), "\n")
cat("Min modal probability:", round(properties$min_max_probability, 3), "\n")
cat("Effective range:", paste(round(properties$effective_range, 2), collapse = " to "), "\n")
# Quality assessment
quality_score <- 0
# All categories used
if(properties$n_categories_used == length(examples[[i]]$thresholds) + 1) {
quality_score <- quality_score + 2
}
# High discrimination
if(examples[[i]]$discrimination > 1.0) {
quality_score <- quality_score + 1
}
# Good modal probabilities
if(properties$min_max_probability > 0.3) {
quality_score <- quality_score + 1
}
# Well-spaced thresholds
threshold_spacings <- diff(examples[[i]]$thresholds)
if(all(threshold_spacings > 0.5) && all(threshold_spacings < 2.5)) {
quality_score <- quality_score + 1
}
quality_assessment <- ifelse(quality_score >= 4, "Excellent",
ifelse(quality_score >= 3, "Good",
ifelse(quality_score >= 2, "Fair", "Poor")))
cat("Quality assessment:", quality_assessment, "(", quality_score, "/5 )\n\n")
}
return(list(
examples = examples,
crf_data = all_crf_data,
theta_range = theta_range
))
}
# Run CRF demonstration
crf_results <- demonstrate_category_response_functions()## === CATEGORY RESPONSE FUNCTIONS DEMONSTRASIYASI ===
## Creating Category Response Function plots...
##
## === CRF CHARACTERISTICS ANALYSIS ===
## Analysis for Well-Functioning Item :
## Categories used: 0, 1, 2, 3, 4
## Crossing points: -1.11, -0.5, 0.5, 1.11
## Min modal probability: 0.318
## Effective range: -4 to 4
## Quality assessment: Excellent ( 5 /5 )
##
## Analysis for Low Discrimination Item :
## Categories used: 0, 4
## Crossing points: -0.5, -0.5, 0.5, 0.5
## Min modal probability: 0.231
## Effective range: -4 to 4
## Quality assessment: Poor ( 1 /5 )
##
## Analysis for High Discrimination Item :
## Categories used: 0, 1, 2, 3, 4
## Crossing points: -1.43, -0.5, 0.5, 1.43
## Min modal probability: 0.424
## Effective range: -4 to 4
## Quality assessment: Excellent ( 5 /5 )
##
## Analysis for Poorly Spaced Thresholds :
## Categories used: 0, 4
## Crossing points: -1, -1.8, -1.6, -2.4
## Min modal probability: 0.389
## Effective range: -4 to 4
## Quality assessment: Fair ( 2 /5 )##
## === CRF INTERPRETATION GUIDELINES ===# CRF interpretation guidelines
crf_guidelines <- c(
"**IDEAL CRF CHARACTERISTICS:**",
"",
"1. **Category Usage:**",
" • All categories should be modal somewhere on θ scale",
" • Each category should have peak probability > 0.3",
" • No 'dead' categories that are never most likely",
"",
"2. **Threshold Spacing:**",
" • Adjacent thresholds spaced 0.5 to 2.0 logits apart",
" • Equal spacing often optimal for rating scales",
" • Avoid clustering of thresholds",
"",
"3. **Discrimination Effects:**",
" • Higher discrimination → steeper CRFs",
" • Sharper transitions between categories",
" • More precise ability discrimination",
"",
"4. **Troubleshooting Poor CRFs:**",
" • Flat CRFs → increase discrimination or collapse categories",
" • Unused categories → adjust thresholds or reduce categories",
" • Disordered thresholds → violates model assumptions"
)
for(guideline in crf_guidelines) {
cat(guideline, "\n")
}## **IDEAL CRF CHARACTERISTICS:**
##
## 1. **Category Usage:**
## • All categories should be modal somewhere on θ scale
## • Each category should have peak probability > 0.3
## • No 'dead' categories that are never most likely
##
## 2. **Threshold Spacing:**
## • Adjacent thresholds spaced 0.5 to 2.0 logits apart
## • Equal spacing often optimal for rating scales
## • Avoid clustering of thresholds
##
## 3. **Discrimination Effects:**
## • Higher discrimination → steeper CRFs
## • Sharper transitions between categories
## • More precise ability discrimination
##
## 4. **Troubleshooting Poor CRFs:**
## • Flat CRFs → increase discrimination or collapse categories
## • Unused categories → adjust thresholds or reduce categories
## • Disordered thresholds → violates model assumptions##
## === CRF QUALITY INDICATORS ===# Quality indicators table
quality_indicators <- data.frame(
Indicator = c("Category Usage", "Peak Probabilities", "Threshold Spacing",
"Discrimination Level", "Effective Range", "Model Fit"),
Good_Range = c("All categories used", "> 0.30", "0.5 - 2.0 logits",
"> 1.0", "Wide coverage", "Adequate"),
Warning_Signs = c("Unused categories", "< 0.20", "< 0.5 or > 2.5",
"< 0.7", "Narrow range", "Poor fit"),
Remedial_Actions = c("Collapse categories", "Adjust thresholds", "Redistribute thresholds",
"Check item quality", "Review content", "Model modification"),
Software_Checks = c("Modal category plots", "CRF maximum values", "Threshold differences",
"Parameter estimates", "Information plots", "Fit statistics")
)
print(quality_indicators)## Indicator Good_Range Warning_Signs
## 1 Category Usage All categories used Unused categories
## 2 Peak Probabilities > 0.30 < 0.20
## 3 Threshold Spacing 0.5 - 2.0 logits < 0.5 or > 2.5
## 4 Discrimination Level > 1.0 < 0.7
## 5 Effective Range Wide coverage Narrow range
## 6 Model Fit Adequate Poor fit
## Remedial_Actions Software_Checks
## 1 Collapse categories Modal category plots
## 2 Adjust thresholds CRF maximum values
## 3 Redistribute thresholds Threshold differences
## 4 Check item quality Parameter estimates
## 5 Review content Information plots
## 6 Model modification Fit statistics4 Threshold Parameters
4.1 Understanding və Interpretation
## === THRESHOLD PARAMETERS ===## Threshold Parameters: GRM-də kateqoriya keçidlərini müəyyən edən parameterlər## • bik: i item-inin k threshold-u## • P(X ≥ k | θ = bik) = 0.5 interpretasiyası## • Threshold ordering və constraints## • Praktik interpretasiya və meaningful differences# Comprehensive threshold parameters demonstration
demonstrate_threshold_parameters <- function() {
cat("=== THRESHOLD PARAMETERS DEMONSTRASIYASI ===\n")
# Create examples with different threshold configurations
theta_range <- seq(-3, 3, 0.1)
# Example scenarios for threshold analysis
threshold_scenarios <- list(
"Well_Spaced" = list(
discrimination = 1.5,
thresholds = c(-1.5, -0.5, 0.5, 1.5),
description = "Well-spaced thresholds"
),
"Close_Thresholds" = list(
discrimination = 1.5,
thresholds = c(-0.3, -0.1, 0.1, 0.3),
description = "Closely spaced thresholds"
),
"Wide_Thresholds" = list(
discrimination = 1.5,
thresholds = c(-2.5, -1.0, 1.0, 2.5),
description = "Widely spaced thresholds"
),
"Unequal_Spacing" = list(
discrimination = 1.5,
thresholds = c(-2.0, -1.8, 0.5, 2.0),
description = "Unequally spaced thresholds"
)
)
cat("Threshold Scenarios Created:\n")
for(name in names(threshold_scenarios)) {
scenario <- threshold_scenarios[[name]]
cat(name, ":", scenario$description, "\n")
cat(" Thresholds:", paste(scenario$thresholds, collapse = ", "), "\n")
}
cat("\n")
# Analyze threshold properties
analyze_threshold_properties <- function(scenario_name, params) {
cat("=== ANALYSIS:", scenario_name, "===\n")
# Basic threshold statistics
thresholds <- params$thresholds
n_thresholds <- length(thresholds)
# Threshold spacing
threshold_spacing <- diff(thresholds)
cat("Threshold Analysis:\n")
cat("Number of thresholds:", n_thresholds, "\n")
cat("Threshold values:", paste(round(thresholds, 2), collapse = ", "), "\n")
cat("Threshold spacings:", paste(round(threshold_spacing, 2), collapse = ", "), "\n")
cat("Mean spacing:", round(mean(threshold_spacing), 2), "\n")
cat("SD spacing:", round(sd(threshold_spacing), 2), "\n")
# Expected score at thresholds
expected_scores_at_thresholds <- numeric(n_thresholds)
for(i in 1:n_thresholds) {
theta_threshold <- thresholds[i]
cat_probs <- grm_probability(theta_threshold, params$discrimination, thresholds)
expected_scores_at_thresholds[i] <- sum(cat_probs * (0:(length(cat_probs) - 1)))
}
cat("Expected scores at thresholds:", paste(round(expected_scores_at_thresholds, 2), collapse = ", "), "\n")
# Calculate cumulative probabilities at thresholds
cumulative_probs_at_thresholds <- matrix(NA, n_thresholds, n_thresholds)
for(i in 1:n_thresholds) {
theta_threshold <- thresholds[i]
for(k in 1:n_thresholds) {
cum_prob <- exp(params$discrimination * (theta_threshold - thresholds[k])) /
(1 + exp(params$discrimination * (theta_threshold - thresholds[k])))
cumulative_probs_at_thresholds[i, k] <- cum_prob
}
}
cat("Cumulative probabilities at own threshold (should be ≈ 0.5):\n")
for(i in 1:n_thresholds) {
cat(" Threshold", i, ":", round(cumulative_probs_at_thresholds[i, i], 3), "\n")
}
# Category probabilities at each threshold
cat("\nCategory probabilities at each threshold:\n")
for(i in 1:n_thresholds) {
theta_threshold <- thresholds[i]
cat_probs <- grm_probability(theta_threshold, params$discrimination, thresholds)
cat("At threshold", i, "(θ =", round(theta_threshold, 2), "):\n")
for(k in 1:length(cat_probs)) {
cat(" P(X =", k-1, ") =", round(cat_probs[k], 3), "\n")
}
}
# Quality assessment
quality_indicators <- list()
# Spacing uniformity
spacing_cv <- sd(threshold_spacing) / mean(threshold_spacing)
quality_indicators$spacing_uniformity <- ifelse(spacing_cv < 0.3, "Good", "Poor")
# Adequate spacing
min_spacing <- min(threshold_spacing)
quality_indicators$adequate_spacing <- ifelse(min_spacing > 0.5, "Good", "Concern")
# Not too wide
max_spacing <- max(threshold_spacing)
quality_indicators$not_too_wide <- ifelse(max_spacing < 2.5, "Good", "Concern")
# Expected score progression
score_increases <- diff(expected_scores_at_thresholds)
quality_indicators$score_progression <- ifelse(all(score_increases > 0), "Good", "Problem")
cat("\nQuality Assessment:\n")
for(indicator in names(quality_indicators)) {
cat(indicator, ":", quality_indicators[[indicator]], "\n")
}
cat("\n")
return(list(
thresholds = thresholds,
spacing = threshold_spacing,
expected_scores = expected_scores_at_thresholds,
quality = quality_indicators
))
}
# Analyze all scenarios
scenario_analyses <- list()
for(scenario_name in names(threshold_scenarios)) {
analysis <- analyze_threshold_properties(scenario_name, threshold_scenarios[[scenario_name]])
scenario_analyses[[scenario_name]] <- analysis
}
# Create threshold visualization
if(require(ggplot2, quietly = TRUE)) {
cat("=== THRESHOLD VISUALIZATION ===\n")
# Create cumulative probability plots
for(scenario_name in names(threshold_scenarios)) {
params <- threshold_scenarios[[scenario_name]]
# Calculate cumulative probabilities
cum_prob_data <- data.frame()
for(theta in theta_range) {
for(k in 1:length(params$thresholds)) {
cum_prob <- exp(params$discrimination * (theta - params$thresholds[k])) /
(1 + exp(params$discrimination * (theta - params$thresholds[k])))
cum_prob_data <- rbind(cum_prob_data, data.frame(
Theta = theta,
Threshold = paste0("P*(≥", k, ")"),
Cumulative_Probability = cum_prob,
Scenario = params$description
))
}
}
# Create plot
p <- ggplot(cum_prob_data, aes(x = Theta, y = Cumulative_Probability, color = Threshold)) +
geom_line(size = 1.2) +
geom_hline(yintercept = 0.5, linetype = "dashed", alpha = 0.5) +
geom_vline(xintercept = params$thresholds, linetype = "dotted", alpha = 0.7) +
labs(title = paste("Cumulative Probabilities:", params$description),
subtitle = paste("Thresholds:", paste(round(params$thresholds, 2), collapse = ", ")),
x = "Ability (θ)",
y = "Cumulative Probability P*(X ≥ k)",
color = "Threshold") +
scale_color_brewer(type = "qual", palette = "Dark2") +
theme_minimal() +
theme(plot.title = element_text(hjust = 0.5, size = 14, face = "bold"),
plot.subtitle = element_text(hjust = 0.5, size = 11),
legend.position = "bottom") +
ylim(0, 1) +
xlim(-3, 3)
print(p)
}
}
# Threshold interpretation guidelines
cat("=== THRESHOLD INTERPRETATION GUIDELINES ===\n")
interpretation_guidelines <- c(
"**THRESHOLD PARAMETER MEANING:**",
"",
"bik = Ability level where P(X ≥ k | θ) = 0.5",
"",
"• Threshold 1 (bi1): 50% chance of scoring 1 or higher",
"• Threshold 2 (bi2): 50% chance of scoring 2 or higher",
"• And so on...",
"",
"**PRACTICAL INTERPRETATION:**",
"",
"• Lower thresholds → easier to achieve higher categories",
"• Higher thresholds → harder to achieve higher categories",
"• Threshold spacing affects category usage patterns",
"",
"**OPTIMAL THRESHOLD SPACING:**",
"",
"• Generally 0.5 to 2.0 logits apart",
"• Equal spacing often works well for rating scales",
"• Closer spacing → more categories used near center",
"• Wider spacing → broader category usage",
"",
"**THRESHOLD PROBLEMS:**",
"",
"• Disordered thresholds (bi1 > bi2): Model violation",
"• Very close thresholds (< 0.3): Categories hard to distinguish",
"• Very distant thresholds (> 3.0): Missing intermediate responses"
)
for(guideline in interpretation_guidelines) {
cat(guideline, "\n")
}
return(list(
scenarios = threshold_scenarios,
analyses = scenario_analyses
))
}
# Run threshold parameters demonstration
threshold_results <- demonstrate_threshold_parameters()## === THRESHOLD PARAMETERS DEMONSTRASIYASI ===
## Threshold Scenarios Created:
## Well_Spaced : Well-spaced thresholds
## Thresholds: -1.5, -0.5, 0.5, 1.5
## Close_Thresholds : Closely spaced thresholds
## Thresholds: -0.3, -0.1, 0.1, 0.3
## Wide_Thresholds : Widely spaced thresholds
## Thresholds: -2.5, -1, 1, 2.5
## Unequal_Spacing : Unequally spaced thresholds
## Thresholds: -2, -1.8, 0.5, 2
##
## === ANALYSIS: Well_Spaced ===
## Threshold Analysis:
## Number of thresholds: 4
## Threshold values: -1.5, -0.5, 0.5, 1.5
## Threshold spacings: 1, 1, 1
## Mean spacing: 1
## SD spacing: 0
## Expected scores at thresholds: 0.74, 1.55, 2.45, 3.26
## Cumulative probabilities at own threshold (should be ≈ 0.5):
## Threshold 1 : 0.5
## Threshold 2 : 0.5
## Threshold 3 : 0.5
## Threshold 4 : 0.5
##
## Category probabilities at each threshold:
## At threshold 1 (θ = -1.5 ):
## P(X = 0 ) = 0.5
## P(X = 1 ) = 0.318
## P(X = 2 ) = 0.135
## P(X = 3 ) = 0.036
## P(X = 4 ) = 0.011
## At threshold 2 (θ = -0.5 ):
## P(X = 0 ) = 0.182
## P(X = 1 ) = 0.318
## P(X = 2 ) = 0.318
## P(X = 3 ) = 0.135
## P(X = 4 ) = 0.047
## At threshold 3 (θ = 0.5 ):
## P(X = 0 ) = 0.047
## P(X = 1 ) = 0.135
## P(X = 2 ) = 0.318
## P(X = 3 ) = 0.318
## P(X = 4 ) = 0.182
## At threshold 4 (θ = 1.5 ):
## P(X = 0 ) = 0.011
## P(X = 1 ) = 0.036
## P(X = 2 ) = 0.135
## P(X = 3 ) = 0.318
## P(X = 4 ) = 0.5
##
## Quality Assessment:
## spacing_uniformity : Good
## adequate_spacing : Good
## not_too_wide : Good
## score_progression : Good
##
## === ANALYSIS: Close_Thresholds ===
## Threshold Analysis:
## Number of thresholds: 4
## Threshold values: -0.3, -0.1, 0.1, 0.3
## Threshold spacings: 0.2, 0.2, 0.2
## Mean spacing: 0.2
## SD spacing: 0
## Expected scores at thresholds: 1.57, 1.85, 2.15, 2.43
## Cumulative probabilities at own threshold (should be ≈ 0.5):
## Threshold 1 : 0.5
## Threshold 2 : 0.5
## Threshold 3 : 0.5
## Threshold 4 : 0.5
##
## Category probabilities at each threshold:
## At threshold 1 (θ = -0.3 ):
## P(X = 0 ) = 0.5
## P(X = 1 ) = 0.074
## P(X = 2 ) = 0.071
## P(X = 3 ) = 0.065
## P(X = 4 ) = 0.289
## At threshold 2 (θ = -0.1 ):
## P(X = 0 ) = 0.426
## P(X = 1 ) = 0.074
## P(X = 2 ) = 0.074
## P(X = 3 ) = 0.071
## P(X = 4 ) = 0.354
## At threshold 3 (θ = 0.1 ):
## P(X = 0 ) = 0.354
## P(X = 1 ) = 0.071
## P(X = 2 ) = 0.074
## P(X = 3 ) = 0.074
## P(X = 4 ) = 0.426
## At threshold 4 (θ = 0.3 ):
## P(X = 0 ) = 0.289
## P(X = 1 ) = 0.065
## P(X = 2 ) = 0.071
## P(X = 3 ) = 0.074
## P(X = 4 ) = 0.5
##
## Quality Assessment:
## spacing_uniformity : Good
## adequate_spacing : Concern
## not_too_wide : Good
## score_progression : Good
##
## === ANALYSIS: Wide_Thresholds ===
## Threshold Analysis:
## Number of thresholds: 4
## Threshold values: -2.5, -1, 1, 2.5
## Threshold spacings: 1.5, 2, 1.5
## Mean spacing: 1.67
## SD spacing: 0.29
## Expected scores at thresholds: 0.6, 1.46, 2.54, 3.4
## Cumulative probabilities at own threshold (should be ≈ 0.5):
## Threshold 1 : 0.5
## Threshold 2 : 0.5
## Threshold 3 : 0.5
## Threshold 4 : 0.5
##
## Category probabilities at each threshold:
## At threshold 1 (θ = -2.5 ):
## P(X = 0 ) = 0.5
## P(X = 1 ) = 0.405
## P(X = 2 ) = 0.09
## P(X = 3 ) = 0.005
## P(X = 4 ) = 0.001
## At threshold 2 (θ = -1 ):
## P(X = 0 ) = 0.095
## P(X = 1 ) = 0.405
## P(X = 2 ) = 0.453
## P(X = 3 ) = 0.042
## P(X = 4 ) = 0.005
## At threshold 3 (θ = 1 ):
## P(X = 0 ) = 0.005
## P(X = 1 ) = 0.042
## P(X = 2 ) = 0.453
## P(X = 3 ) = 0.405
## P(X = 4 ) = 0.095
## At threshold 4 (θ = 2.5 ):
## P(X = 0 ) = 0.001
## P(X = 1 ) = 0.005
## P(X = 2 ) = 0.09
## P(X = 3 ) = 0.405
## P(X = 4 ) = 0.5
##
## Quality Assessment:
## spacing_uniformity : Good
## adequate_spacing : Good
## not_too_wide : Good
## score_progression : Good
##
## === ANALYSIS: Unequal_Spacing ===
## Threshold Analysis:
## Number of thresholds: 4
## Threshold values: -2, -1.8, 0.5, 2
## Threshold spacings: 0.2, 2.3, 1.5
## Mean spacing: 1.33
## SD spacing: 1.06
## Expected scores at thresholds: 0.95, 1.11, 2.54, 3.4
## Cumulative probabilities at own threshold (should be ≈ 0.5):
## Threshold 1 : 0.5
## Threshold 2 : 0.5
## Threshold 3 : 0.5
## Threshold 4 : 0.5
##
## Category probabilities at each threshold:
## At threshold 1 (θ = -2 ):
## P(X = 0 ) = 0.5
## P(X = 1 ) = 0.074
## P(X = 2 ) = 0.403
## P(X = 3 ) = 0.021
## P(X = 4 ) = 0.002
## At threshold 2 (θ = -1.8 ):
## P(X = 0 ) = 0.426
## P(X = 1 ) = 0.074
## P(X = 2 ) = 0.469
## P(X = 3 ) = 0.027
## P(X = 4 ) = 0.003
## At threshold 3 (θ = 0.5 ):
## P(X = 0 ) = 0.023
## P(X = 1 ) = 0.008
## P(X = 2 ) = 0.469
## P(X = 3 ) = 0.405
## P(X = 4 ) = 0.095
## At threshold 4 (θ = 2 ):
## P(X = 0 ) = 0.002
## P(X = 1 ) = 0.001
## P(X = 2 ) = 0.092
## P(X = 3 ) = 0.405
## P(X = 4 ) = 0.5
##
## Quality Assessment:
## spacing_uniformity : Poor
## adequate_spacing : Concern
## not_too_wide : Good
## score_progression : Good
##
## === THRESHOLD VISUALIZATION ===
## === THRESHOLD INTERPRETATION GUIDELINES ===
## **THRESHOLD PARAMETER MEANING:**
##
## bik = Ability level where P(X ≥ k | θ) = 0.5
##
## • Threshold 1 (bi1): 50% chance of scoring 1 or higher
## • Threshold 2 (bi2): 50% chance of scoring 2 or higher
## • And so on...
##
## **PRACTICAL INTERPRETATION:**
##
## • Lower thresholds → easier to achieve higher categories
## • Higher thresholds → harder to achieve higher categories
## • Threshold spacing affects category usage patterns
##
## **OPTIMAL THRESHOLD SPACING:**
##
## • Generally 0.5 to 2.0 logits apart
## • Equal spacing often works well for rating scales
## • Closer spacing → more categories used near center
## • Wider spacing → broader category usage
##
## **THRESHOLD PROBLEMS:**
##
## • Disordered thresholds (bi1 > bi2): Model violation
## • Very close thresholds (< 0.3): Categories hard to distinguish
## • Very distant thresholds (> 3.0): Missing intermediate responses##
## === THRESHOLD PARAMETER ESTIMATION ===# Parameter estimation considerations
estimation_considerations <- data.frame(
Aspect = c("Sample Size", "Category Frequency", "Threshold Stability",
"Convergence", "Starting Values", "Constraints"),
Requirement = c("50+ per category", "5+ responses per category", "Large N for precision",
"Proper algorithm", "Reasonable initial values", "Ordered thresholds"),
Common_Problems = c("Small samples", "Empty categories", "Boundary estimates",
"Non-convergence", "Poor starting values", "Disordered estimates"),
Solutions = c("Increase N or collapse", "Combine categories", "Add constraints",
"Change algorithm", "Use data-based starts", "Fix ordering"),
Software_Notes = c("mirt, ltm handle automatically", "Automatic collapsing", "Bayesian priors help",
"Multiple optimizers", "Data-driven defaults", "Built-in ordering")
)
print(estimation_considerations)## Aspect Requirement Common_Problems
## 1 Sample Size 50+ per category Small samples
## 2 Category Frequency 5+ responses per category Empty categories
## 3 Threshold Stability Large N for precision Boundary estimates
## 4 Convergence Proper algorithm Non-convergence
## 5 Starting Values Reasonable initial values Poor starting values
## 6 Constraints Ordered thresholds Disordered estimates
## Solutions Software_Notes
## 1 Increase N or collapse mirt, ltm handle automatically
## 2 Combine categories Automatic collapsing
## 3 Add constraints Bayesian priors help
## 4 Change algorithm Multiple optimizers
## 5 Use data-based starts Data-driven defaults
## 6 Fix ordering Built-in ordering##
## === THRESHOLD REPORTING STANDARDS ===# Reporting standards
reporting_standards <- c(
"**REQUIRED REPORTING:**",
"",
"1. **Parameter Estimates:**",
" • All threshold estimates with standard errors",
" • Discrimination parameters",
" • Model identification constraints used",
"",
"2. **Diagnostic Information:**",
" • Threshold spacing analysis",
" • Category usage frequencies",
" • Parameter stability across samples",
"",
"3. **Quality Indicators:**",
" • Model fit statistics",
" • Category response function plots",
" • Information function curves",
"",
"**INTERPRETATION GUIDELINES:**",
"• Relate thresholds to practical meaning",
"• Discuss category spacing implications",
"• Address any problematic patterns",
"• Provide guidance for score interpretation"
)
for(standard in reporting_standards) {
cat(standard, "\n")
}## **REQUIRED REPORTING:**
##
## 1. **Parameter Estimates:**
## • All threshold estimates with standard errors
## • Discrimination parameters
## • Model identification constraints used
##
## 2. **Diagnostic Information:**
## • Threshold spacing analysis
## • Category usage frequencies
## • Parameter stability across samples
##
## 3. **Quality Indicators:**
## • Model fit statistics
## • Category response function plots
## • Information function curves
##
## **INTERPRETATION GUIDELINES:**
## • Relate thresholds to practical meaning
## • Discuss category spacing implications
## • Address any problematic patterns
## • Provide guidance for score interpretation5 Rating Scale Model və Partial Credit Model
5.1 Alternative Polytomous IRT Models
## === RATING SCALE MODEL VƏ PARTIAL CREDIT MODEL ===## Alternative Polytomous Models: GRM-ə alternativ yanaşmalar## • Rating Scale Model (RSM): Ümumi threshold structure## • Partial Credit Model (PCM): Item-specific thresholds## • Generalized Partial Credit Model (GPCM): PCM + discrimination## • Model comparison və selection criteria# Comprehensive RSM and PCM demonstration
demonstrate_rsm_pcm_models <- function() {
cat("=== RSM VƏ PCM MODELS DEMONSTRASIYASI ===\n")
# Generate data suitable for different models
set.seed(1801)
n_items <- 8
n_persons <- 1200
n_categories <- 4 # 0, 1, 2, 3
# Generate person abilities
theta <- rnorm(n_persons, 0, 1)
# Generate item difficulties (location parameters)
item_difficulties <- rnorm(n_items, 0, 1)
# Model 1: Rating Scale Model (RSM)
# Common threshold structure across items
common_thresholds <- c(-1.2, 0.0, 1.2) # τ1, τ2, τ3
# Model 2: Partial Credit Model (PCM)
# Item-specific thresholds
item_thresholds <- matrix(NA, n_items, n_categories - 1)
for(i in 1:n_items) {
# Generate ordered thresholds for each item
item_thresholds[i, ] <- sort(rnorm(n_categories - 1, 0, 0.8))
}
cat("Model Configurations:\n")
cat("RSM - Common thresholds:", paste(round(common_thresholds, 2), collapse = ", "), "\n")
cat("PCM - Item-specific thresholds (first 3 items):\n")
for(i in 1:3) {
cat(" Item", i, ":", paste(round(item_thresholds[i, ], 2), collapse = ", "), "\n")
}
cat("\n")
# RSM probability function
rsm_probability <- function(theta, item_difficulty, common_thresholds) {
n_cats <- length(common_thresholds) + 1
# Calculate numerators for each category
numerators <- numeric(n_cats)
for(k in 1:n_cats) {
if(k == 1) {
numerators[k] <- 1 # Category 0
} else {
sum_thresholds <- sum(common_thresholds[1:(k-1)])
numerators[k] <- exp((k-1) * theta - (k-1) * item_difficulty - sum_thresholds)
}
}
# Denominator
denominator <- sum(numerators)
# Category probabilities
probabilities <- numerators / denominator
return(probabilities)
}
# PCM probability function
pcm_probability <- function(theta, item_difficulty, item_thresholds) {
n_cats <- length(item_thresholds) + 1
# Calculate numerators for each category
numerators <- numeric(n_cats)
for(k in 1:n_cats) {
if(k == 1) {
numerators[k] <- 1 # Category 0
} else {
sum_thresholds <- sum(item_thresholds[1:(k-1)])
numerators[k] <- exp((k-1) * theta - (k-1) * item_difficulty - sum_thresholds)
}
}
# Denominator
denominator <- sum(numerators)
# Category probabilities
probabilities <- numerators / denominator
return(probabilities)
}
# Generate data under RSM
rsm_responses <- matrix(NA, n_persons, n_items)
for(i in 1:n_persons) {
for(j in 1:n_items) {
cat_probs <- rsm_probability(theta[i], item_difficulties[j], common_thresholds)
rsm_responses[i, j] <- sample(0:(n_categories-1), 1, prob = cat_probs)
}
}
# Generate data under PCM
pcm_responses <- matrix(NA, n_persons, n_items)
for(i in 1:n_persons) {
for(j in 1:n_items) {
cat_probs <- pcm_probability(theta[i], item_difficulties[j], item_thresholds[j, ])
pcm_responses[i, j] <- sample(0:(n_categories-1), 1, prob = cat_probs)
}
}
colnames(rsm_responses) <- paste0("Item", 1:n_items)
colnames(pcm_responses) <- paste0("Item", 1:n_items)
cat("Response Data Generated:\n")
cat("RSM data range:", range(rsm_responses), "\n")
cat("PCM data range:", range(pcm_responses), "\n")
# Descriptive analysis
cat("\n=== DESCRIPTIVE ANALYSIS ===\n")
# Category usage patterns
analyze_category_usage <- function(responses, model_name) {
cat("Category Usage -", model_name, ":\n")
overall_usage <- table(factor(as.vector(responses), levels = 0:(n_categories-1)))
cat("Overall:", paste(overall_usage, collapse = " / "),
"(", paste(round(prop.table(overall_usage) * 100, 1), collapse = "% / "), "%)\n")
# Item-specific usage
item_usage <- data.frame()
for(j in 1:n_items) {
item_cats <- table(factor(responses[, j], levels = 0:(n_categories-1)))
item_usage <- rbind(item_usage, data.frame(
Item = j,
Cat0 = item_cats[1],
Cat1 = item_cats[2],
Cat2 = item_cats[3],
Cat3 = item_cats[4],
Mean = round(mean(responses[, j]), 2),
SD = round(sd(responses[, j]), 2)
))
}
cat("Item-level usage (first 4 items):\n")
print(item_usage[1:4, ])
return(item_usage)
}
rsm_usage <- analyze_category_usage(rsm_responses, "RSM")
cat("\n")
pcm_usage <- analyze_category_usage(pcm_responses, "PCM")
# Model comparison using available packages
cat("\n=== MODEL FITTING VƏ COMPARISON ===\n")
if(require(eRm, quietly = TRUE)) {
cat("Fitting models with eRm package...\n")
tryCatch({
# Fit RSM
rsm_fit <- RSM(rsm_responses)
cat("RSM fitted successfully\n")
# Fit PCM
pcm_fit <- PCM(pcm_responses)
cat("PCM fitted successfully\n")
# Model comparison would go here
# Note: eRm uses different parameterization than our simulation
}, error = function(e) {
cat("eRm fitting failed:", e$message, "\n")
})
}
if(require(mirt, quietly = TRUE)) {
cat("\nFitting models with mirt package...\n")
tryCatch({
# Fit RSM (constrained GPCM)
rsm_mirt <- mirt(rsm_responses, 1, itemtype = "rsm", verbose = FALSE)
cat("RSM fitted with mirt\n")
# Fit PCM (partial credit model)
pcm_mirt <- mirt(pcm_responses, 1, itemtype = "Rasch", verbose = FALSE)
cat("PCM fitted with mirt\n")
# Fit GPCM (generalized partial credit)
gpcm_mirt <- mirt(pcm_responses, 1, itemtype = "gpcm", verbose = FALSE)
cat("GPCM fitted with mirt\n")
# Extract fit statistics
rsm_fit_stats <- data.frame(
Model = "RSM",
AIC = extract.mirt(rsm_mirt, "AIC"),
BIC = extract.mirt(rsm_mirt, "BIC"),
LogLik = extract.mirt(rsm_mirt, "logLik")
)
pcm_fit_stats <- data.frame(
Model = "PCM",
AIC = extract.mirt(pcm_mirt, "AIC"),
BIC = extract.mirt(pcm_mirt, "BIC"),
LogLik = extract.mirt(pcm_mirt, "logLik")
)
gpcm_fit_stats <- data.frame(
Model = "GPCM",
AIC = extract.mirt(gpcm_mirt, "AIC"),
BIC = extract.mirt(gpcm_mirt, "BIC"),
LogLik = extract.mirt(gpcm_mirt, "logLik")
)
fit_comparison <- rbind(rsm_fit_stats, pcm_fit_stats, gpcm_fit_stats)
fit_comparison$Delta_AIC <- fit_comparison$AIC - min(fit_comparison$AIC)
fit_comparison$Delta_BIC <- fit_comparison$BIC - min(fit_comparison$BIC)
cat("\nModel Fit Comparison:\n")
print(fit_comparison)
# Parameter comparison
cat("\nParameter Comparison (first 3 items):\n")
rsm_coefs <- coef(rsm_mirt, simplify = TRUE)
pcm_coefs <- coef(pcm_mirt, simplify = TRUE)
gpcm_coefs <- coef(gpcm_mirt, simplify = TRUE)
cat("RSM item parameters (sample):\n")
print(round(rsm_coefs$items[1:3, ], 3))
cat("\nGPCM item parameters (sample):\n")
print(round(gpcm_coefs$items[1:3, ], 3))
}, error = function(e) {
cat("mirt fitting failed:", e$message, "\n")
fit_comparison <- NULL
})
}
# Model characteristic comparison
cat("\n=== MODEL CHARACTERISTICS COMPARISON ===\n")
model_comparison <- data.frame(
Aspect = c("Parameterization", "Threshold Structure", "Discrimination",
"Complexity", "Assumptions", "Best For"),
GRM = c("Cumulative logits", "Item-specific ordered", "Item-specific",
"Most complex", "Unidimensional", "Likert scales"),
RSM = c("Adjacent categories", "Common across items", "Equal (=1)",
"Least complex", "Rasch family", "Rating scales"),
PCM = c("Adjacent categories", "Item-specific", "Equal (=1)",
"Moderate", "Rasch family", "Achievement tests"),
GPCM = c("Adjacent categories", "Item-specific", "Item-specific",
"Complex", "Unidimensional", "Mixed item types")
)
print(model_comparison)
return(list(
rsm_data = rsm_responses,
pcm_data = pcm_responses,
true_parameters = list(
theta = theta,
item_difficulties = item_difficulties,
common_thresholds = common_thresholds,
item_thresholds = item_thresholds
),
usage_analysis = list(rsm = rsm_usage, pcm = pcm_usage),
fit_comparison = if(exists("fit_comparison")) fit_comparison else NULL
))
}
# Run RSM and PCM demonstration
rsm_pcm_results <- demonstrate_rsm_pcm_models()## === RSM VƏ PCM MODELS DEMONSTRASIYASI ===
## Model Configurations:
## RSM - Common thresholds: -1.2, 0, 1.2
## PCM - Item-specific thresholds (first 3 items):
## Item 1 : -0.29, 0.58, 0.66
## Item 2 : -0.82, -0.59, -0.1
## Item 3 : -0.56, -0.11, 1.54
##
## Response Data Generated:
## RSM data range: 0 3
## PCM data range: 0 3
##
## === DESCRIPTIVE ANALYSIS ===
## Category Usage - RSM :
## Overall: 1940 / 2535 / 2723 / 2402 ( 20.2% / 26.4% / 28.4% / 25 %)
## Item-level usage (first 4 items):
## Item Cat0 Cat1 Cat2 Cat3 Mean SD
## 0 1 25 136 387 652 2.39 0.77
## 01 2 276 410 346 168 1.34 0.98
## 02 3 235 417 350 198 1.43 0.98
## 03 4 103 299 434 364 1.88 0.94
##
## Category Usage - PCM :
## Overall: 2368 / 2193 / 2138 / 2901 ( 24.7% / 22.8% / 22.3% / 30.2 %)
## Item-level usage (first 4 items):
## Item Cat0 Cat1 Cat2 Cat3 Mean SD
## 0 1 50 144 256 750 2.42 0.86
## 01 2 279 260 274 387 1.64 1.16
## 02 3 377 306 365 152 1.24 1.03
## 03 4 175 404 306 315 1.63 1.02
##
## === MODEL FITTING VƏ COMPARISON ===
## Fitting models with eRm package...
## RSM fitted successfully
## PCM fitted successfully
##
## Fitting models with mirt package...
## RSM fitted with mirt
## PCM fitted with mirt
## GPCM fitted with mirt
##
## Model Fit Comparison:
## Model AIC BIC LogLik Delta_AIC Delta_BIC
## 1 RSM 20980.31 21036.31 -10479.16 0.0000 0.0000
## 2 PCM 21504.92 21632.17 -10727.46 524.6035 595.8645
## 3 GPCM 21514.12 21677.00 -10725.06 533.8038 640.6955
##
## Parameter Comparison (first 3 items):
## RSM item parameters (sample):
## a1 b1 b2 b3 c
## Item1 1 -3.008 -1.719 -0.482 0.000
## Item2 1 -3.008 -1.719 -0.482 -2.018
## Item3 1 -3.008 -1.719 -0.482 -1.866
##
## GPCM item parameters (sample):
## a1 ak0 ak1 ak2 ak3 d0 d1 d2 d3
## Item1 0.958 0 1 2 3 0 2.201 3.361 4.407
## Item2 0.936 0 1 2 3 0 0.510 0.634 0.545
## Item3 0.992 0 1 2 3 0 0.227 0.242 -1.374
##
## === MODEL CHARACTERISTICS COMPARISON ===
## Aspect GRM RSM
## 1 Parameterization Cumulative logits Adjacent categories
## 2 Threshold Structure Item-specific ordered Common across items
## 3 Discrimination Item-specific Equal (=1)
## 4 Complexity Most complex Least complex
## 5 Assumptions Unidimensional Rasch family
## 6 Best For Likert scales Rating scales
## PCM GPCM
## 1 Adjacent categories Adjacent categories
## 2 Item-specific Item-specific
## 3 Equal (=1) Item-specific
## 4 Moderate Complex
## 5 Rasch family Unidimensional
## 6 Achievement tests Mixed item types##
## === MODEL SELECTION GUIDELINES ===# Model selection guidelines
selection_guidelines <- c(
"**MODEL SELECTION CRITERIA:**",
"",
"1. **Theoretical Considerations:**",
" • GRM: For Likert-type scales with ordered responses",
" • RSM: When threshold structure should be common across items",
" • PCM: For achievement tests with item-specific difficulty progression",
" • GPCM: When items differ in discrimination",
"",
"2. **Data Characteristics:**",
" • Sample size: RSM/PCM need less data than GRM/GPCM",
" • Category usage: All models require adequate category frequencies",
" • Item heterogeneity: More varied items favor GRM/GPCM",
"",
"3. **Statistical Criteria:**",
" • Model fit: Chi-square, AIC, BIC comparisons",
" • Parameter stability: Reasonable estimates and standard errors",
" • Residual analysis: Person and item fit statistics",
"",
"4. **Practical Considerations:**",
" • Interpretability: Simpler models often preferred",
" • Software availability: Not all packages support all models",
" • Reporting requirements: Some contexts prefer specific approaches"
)
for(guideline in selection_guidelines) {
cat(guideline, "\n")
}## **MODEL SELECTION CRITERIA:**
##
## 1. **Theoretical Considerations:**
## • GRM: For Likert-type scales with ordered responses
## • RSM: When threshold structure should be common across items
## • PCM: For achievement tests with item-specific difficulty progression
## • GPCM: When items differ in discrimination
##
## 2. **Data Characteristics:**
## • Sample size: RSM/PCM need less data than GRM/GPCM
## • Category usage: All models require adequate category frequencies
## • Item heterogeneity: More varied items favor GRM/GPCM
##
## 3. **Statistical Criteria:**
## • Model fit: Chi-square, AIC, BIC comparisons
## • Parameter stability: Reasonable estimates and standard errors
## • Residual analysis: Person and item fit statistics
##
## 4. **Practical Considerations:**
## • Interpretability: Simpler models often preferred
## • Software availability: Not all packages support all models
## • Reporting requirements: Some contexts prefer specific approaches##
## === POLYTOMOUS MODEL RECOMMENDATIONS ===# Recommendations table
model_recommendations <- data.frame(
Application = c("Likert Scales", "Rating Scales", "Achievement Tests",
"Mixed Item Types", "Attitude Surveys", "Performance Assessments"),
Primary_Choice = c("GRM", "RSM", "PCM/GPCM", "GPCM", "GRM", "GPCM"),
Alternative = c("GPCM", "GRM", "RSM", "GRM", "RSM", "GRM"),
Key_Consideration = c("Response process", "Common structure", "Item difficulty",
"Item variation", "Response style", "Scoring method"),
Sample_Size = c("Large", "Moderate", "Moderate", "Large", "Moderate", "Large"),
Typical_Categories = c("3-7", "3-5", "2-4", "Variable", "3-7", "3-6")
)
print(model_recommendations)## Application Primary_Choice Alternative Key_Consideration
## 1 Likert Scales GRM GPCM Response process
## 2 Rating Scales RSM GRM Common structure
## 3 Achievement Tests PCM/GPCM RSM Item difficulty
## 4 Mixed Item Types GPCM GRM Item variation
## 5 Attitude Surveys GRM RSM Response style
## 6 Performance Assessments GPCM GRM Scoring method
## Sample_Size Typical_Categories
## 1 Large 3-7
## 2 Moderate 3-5
## 3 Moderate 2-4
## 4 Large Variable
## 5 Moderate 3-7
## 6 Large 3-66 Model Estimation və Practical Applications
6.1 GRM Implementation və Real Data Analysis
## === MODEL ESTIMATION VƏ PRACTICAL APPLICATIONS ===## Practical GRM Implementation: Real data ilə GRM tətbiqi## • Parameter estimation strategies## • Model fit assessment## • Diagnostic procedures## • Score interpretation və reporting# Comprehensive practical applications demonstration
demonstrate_practical_applications <- function() {
cat("=== PRACTICAL APPLICATIONS DEMONSTRASIYASI ===\n")
# Create realistic survey data
set.seed(1901)
# Simulate a psychological scale (e.g., life satisfaction)
n_items <- 12
n_respondents <- 800
n_categories <- 7 # 1-7 Likert scale
# Generate respondent abilities (latent trait: life satisfaction)
life_satisfaction <- rnorm(n_respondents, 0, 1)
# Create items with different characteristics
item_params <- data.frame(
item = 1:n_items,
discrimination = c(1.2, 1.8, 2.1, 1.5, 0.9, 1.7, 2.3, 1.4, 1.9, 1.1, 1.6, 2.0),
content_area = rep(c("Overall", "Family", "Work", "Health"), each = 3)
)
# Generate thresholds (transform to 1-7 scale thinking)
true_thresholds <- matrix(NA, n_items, n_categories - 1)
for(i in 1:n_items) {
# Generate ordered thresholds, then adjust for different response tendencies
base_thresholds <- seq(-2.5, 2.5, length.out = n_categories - 1)
# Add some variation
noise <- rnorm(n_categories - 1, 0, 0.2)
item_thresholds <- sort(base_thresholds + noise)
true_thresholds[i, ] <- item_thresholds
}
# Generate responses (1-7 scale)
responses <- matrix(NA, n_respondents, n_items)
for(i in 1:n_respondents) {
for(j in 1:n_items) {
# Convert to 0-6 for calculation, then add 1
cat_probs <- grm_probability(life_satisfaction[i],
item_params$discrimination[j],
true_thresholds[j, ])
response_0to6 <- sample(0:(n_categories-1), 1, prob = cat_probs)
responses[i, j] <- response_0to6 + 1 # Convert to 1-7 scale
}
}
colnames(responses) <- paste0("LS", sprintf("%02d", 1:n_items))
cat("Life Satisfaction Scale Simulation:\n")
cat("Respondents:", n_respondents, "\n")
cat("Items:", n_items, "\n")
cat("Scale range:", range(responses), "\n")
cat("Content areas:", paste(unique(item_params$content_area), collapse = ", "), "\n\n")
# Descriptive analysis
cat("=== DESCRIPTIVE ANALYSIS ===\n")
# Item descriptives
item_descriptives <- data.frame(
Item = colnames(responses),
Mean = round(apply(responses, 2, mean), 2),
SD = round(apply(responses, 2, sd), 2),
Min = apply(responses, 2, min),
Max = apply(responses, 2, max),
Skewness = round(apply(responses, 2, function(x) {
n <- length(x)
mean_x <- mean(x)
sd_x <- sd(x)
skew <- (n / ((n - 1) * (n - 2))) * sum(((x - mean_x) / sd_x)^3)
return(skew)
}), 2)
)
cat("Item Descriptives (first 6 items):\n")
print(item_descriptives[1:6, ])
# Category frequencies
cat("\nCategory Usage Analysis:\n")
overall_freq <- table(factor(as.vector(responses), levels = 1:7))
overall_prop <- prop.table(overall_freq)
cat("Overall category usage:\n")
for(i in 1:7) {
cat("Category", i, ":", overall_freq[i],
"(", round(overall_prop[i] * 100, 1), "%)\n")
}
# Reliability analysis
cat("\n=== RELIABILITY ANALYSIS ===\n")
if(require(psych, quietly = TRUE)) {
tryCatch({
reliability_results <- alpha(responses)
cat("Cronbach's Alpha:", round(reliability_results$total$raw_alpha, 3), "\n")
cat("Standardized Alpha:", round(reliability_results$total$std.alpha, 3), "\n")
# Item-total correlations
item_total_corr <- cor(responses, rowSums(responses))
cat("Item-total correlations (first 6):",
paste(round(item_total_corr[1:6], 3), collapse = ", "), "\n")
}, error = function(e) {
cat("Reliability analysis failed:", e$message, "\n")
})
}
# Fit GRM using mirt
cat("\n=== GRM MODEL FITTING ===\n")
if(require(mirt, quietly = TRUE)) {
tryCatch({
# Convert to 0-6 scale for mirt (mirt expects 0-based categories)
responses_0to6 <- responses - 1
# Fit GRM
cat("Fitting Graded Response Model...\n")
grm_fit <- mirt(responses_0to6, 1, itemtype = "graded", verbose = FALSE)
cat("GRM fitted successfully\n")
# Extract fit statistics
fit_stats <- data.frame(
Statistic = c("Log-likelihood", "AIC", "BIC", "Number of Parameters"),
Value = c(
round(extract.mirt(grm_fit, "logLik"), 2),
round(extract.mirt(grm_fit, "AIC"), 2),
round(extract.mirt(grm_fit, "BIC"), 2),
extract.mirt(grm_fit, "nfact") * n_items + sum(apply(responses_0to6, 2, max))
)
)
cat("\nModel Fit Statistics:\n")
print(fit_stats)
# Extract parameter estimates
grm_coefs <- coef(grm_fit, simplify = TRUE)
# Item parameters
item_estimates <- data.frame(
Item = colnames(responses),
Discrimination = round(grm_coefs$items[, "a1"], 3),
Threshold_1 = round(grm_coefs$items[, "d1"] / (-grm_coefs$items[, "a1"]), 3),
Threshold_2 = round(grm_coefs$items[, "d2"] / (-grm_coefs$items[, "a1"]), 3),
Threshold_3 = round(grm_coefs$items[, "d3"] / (-grm_coefs$items[, "a1"]), 3),
Content = item_params$content_area
)
cat("\nItem Parameter Estimates (first 6 items):\n")
print(item_estimates[1:6, ])
# Parameter recovery analysis
cat("\n=== PARAMETER RECOVERY ANALYSIS ===\n")
estimated_disc <- grm_coefs$items[, "a1"]
true_disc <- item_params$discrimination
disc_correlation <- cor(estimated_disc, true_disc)
disc_rmse <- sqrt(mean((estimated_disc - true_disc)^2))
cat("Discrimination parameter recovery:\n")
cat("Correlation:", round(disc_correlation, 3), "\n")
cat("RMSE:", round(disc_rmse, 3), "\n")
# Threshold recovery (first threshold only for simplicity)
estimated_thresh1 <- grm_coefs$items[, "d1"] / (-grm_coefs$items[, "a1"])
true_thresh1 <- true_thresholds[, 1]
thresh_correlation <- cor(estimated_thresh1, true_thresh1)
thresh_rmse <- sqrt(mean((estimated_thresh1 - true_thresh1)^2))
cat("First threshold recovery:\n")
cat("Correlation:", round(thresh_correlation, 3), "\n")
cat("RMSE:", round(thresh_rmse, 3), "\n")
# Person parameter estimates (ability scores)
cat("\n=== PERSON PARAMETER ESTIMATION ===\n")
person_scores <- fscores(grm_fit, method = "EAP")
score_correlation <- cor(person_scores[, 1], life_satisfaction)
score_rmse <- sqrt(mean((person_scores[, 1] - life_satisfaction)^2))
cat("Person ability estimation:\n")
cat("Correlation with true θ:", round(score_correlation, 3), "\n")
cat("RMSE:", round(score_rmse, 3), "\n")
# Score distribution
cat("Score distribution:\n")
cat("Mean:", round(mean(person_scores[, 1]), 3), "\n")
cat("SD:", round(sd(person_scores[, 1]), 3), "\n")
cat("Range:", paste(round(range(person_scores[, 1]), 2), collapse = " to "), "\n")
# Model diagnostics
cat("\n=== MODEL DIAGNOSTICS ===\n")
# Item fit
item_fit <- itemfit(grm_fit, method = "S_X2")
if(!is.null(item_fit)) {
significant_misfit <- sum(item_fit$p.S_X2 < 0.01, na.rm = TRUE)
cat("Items with significant misfit (p < 0.01):", significant_misfit, "out of", n_items, "\n")
if(significant_misfit > 0) {
misfit_items <- which(item_fit$p.S_X2 < 0.01)
cat("Misfitting items:", paste(colnames(responses)[misfit_items], collapse = ", "), "\n")
}
}
# Person fit
person_fit <- personfit(grm_fit)
if(!is.null(person_fit)) {
significant_person_misfit <- sum(person_fit$p.Zh < 0.01, na.rm = TRUE)
cat("Persons with significant misfit (p < 0.01):", significant_person_misfit,
"out of", n_respondents, "\n")
}
# Create practical interpretation
cat("\n=== PRACTICAL INTERPRETATION ===\n")
# Score conversion table
theta_levels <- seq(-2, 2, 0.5)
expected_scores <- numeric(length(theta_levels))
for(i in 1:length(theta_levels)) {
item_expected_scores <- numeric(n_items)
for(j in 1:n_items) {
cat_probs <- grm_probability(theta_levels[i],
estimated_disc[j],
grm_coefs$items[j, paste0("d", 1:6)] / (-estimated_disc[j]))
item_expected_scores[j] <- sum(cat_probs * (0:6))
}
expected_scores[i] <- sum(item_expected_scores)
}
score_conversion <- data.frame(
Theta_Level = theta_levels,
Expected_Total_Score = round(expected_scores, 1),
Interpretation = c("Very Low", "Low", "Below Average", "Average",
"Above Average", "High", "Very High", "Extremely High", "Maximum")[1:length(theta_levels)]
)
cat("Score Interpretation Table:\n")
print(score_conversion)
return(list(
responses = responses,
grm_fit = grm_fit,
item_estimates = item_estimates,
person_scores = person_scores,
fit_stats = fit_stats,
recovery = list(
disc_corr = disc_correlation,
thresh_corr = thresh_correlation,
score_corr = score_correlation
),
score_conversion = score_conversion
))
}, error = function(e) {
cat("GRM fitting failed:", e$message, "\n")
return(NULL)
})
} else {
cat("mirt package not available\n")
return(NULL)
}
}
# Run practical applications demonstration
practical_results <- demonstrate_practical_applications()## === PRACTICAL APPLICATIONS DEMONSTRASIYASI ===
## Life Satisfaction Scale Simulation:
## Respondents: 800
## Items: 12
## Scale range: 1 7
## Content areas: Overall, Family, Work, Health
##
## === DESCRIPTIVE ANALYSIS ===
## Item Descriptives (first 6 items):
## Item Mean SD Min Max Skewness
## LS01 LS01 3.90 1.58 1 7 0.07
## LS02 LS02 3.96 1.45 1 7 -0.07
## LS03 LS03 3.98 1.27 1 7 0.22
## LS04 LS04 3.98 1.42 1 7 0.09
## LS05 LS05 4.00 1.84 1 7 0.00
## LS06 LS06 3.82 1.34 1 7 0.06
##
## Category Usage Analysis:
## Overall category usage:
## Category 1 : 601 ( 6.3 %)
## Category 2 : 860 ( 9 %)
## Category 3 : 1998 ( 20.8 %)
## Category 4 : 2781 ( 29 %)
## Category 5 : 1994 ( 20.8 %)
## Category 6 : 865 ( 9 %)
## Category 7 : 501 ( 5.2 %)
##
## === RELIABILITY ANALYSIS ===
## Cronbach's Alpha: 0.883
## Standardized Alpha: 0.89
## Item-total correlations (first 6): 0.594, 0.715, 0.743, 0.667, 0.555, 0.669
##
## === GRM MODEL FITTING ===
## Fitting Graded Response Model...
## GRM fitted successfully
##
## Model Fit Statistics:
## Statistic Value
## 1 Log-likelihood -14688.47
## 2 AIC 29544.94
## 3 BIC 29938.45
## 4 Number of Parameters 84.00
##
## Item Parameter Estimates (first 6 items):
## Item Discrimination Threshold_1 Threshold_2 Threshold_3 Content
## LS01 LS01 1.155 -2.346 -1.740 -0.600 Overall
## LS02 LS02 1.757 -2.086 -1.729 -0.368 Overall
## LS03 LS03 2.151 -2.411 -1.802 -0.451 Overall
## LS04 LS04 1.569 -2.581 -1.636 -0.415 Family
## LS05 LS05 0.927 -2.516 -1.228 -0.878 Family
## LS06 LS06 1.612 -2.381 -1.708 -0.289 Family
##
## === PARAMETER RECOVERY ANALYSIS ===
## Discrimination parameter recovery:
## Correlation: 0.976
## RMSE: 0.091
## First threshold recovery:
## Correlation: 0.754
## RMSE: 0.167
##
## === PERSON PARAMETER ESTIMATION ===
## Person ability estimation:
## Correlation with true θ: 0.952
## RMSE: 0.305
## Score distribution:
## Mean: 0
## SD: 0.952
## Range: -2.82 to 2.51
##
## === MODEL DIAGNOSTICS ===
## GRM fitting failed: unused argument (method = "S_X2")##
## === PRACTICAL IMPLEMENTATION CHECKLIST ===# Implementation checklist
implementation_checklist <- c(
"**PRE-ANALYSIS CHECKS:**",
"☐ Adequate sample size (200+ respondents minimum)",
"☐ Sufficient category usage (5+ responses per category)",
"☐ Data quality (missing data patterns, outliers)",
"☐ Assumption verification (unidimensionality, local independence)",
"",
"**MODEL FITTING:**",
"☐ Software selection and parameterization",
"☐ Starting value specification",
"☐ Convergence monitoring",
"☐ Parameter estimate reasonableness",
"",
"**MODEL EVALUATION:**",
"☐ Overall model fit assessment",
"☐ Item-level fit diagnostics",
"☐ Person-level fit analysis",
"☐ Parameter stability across subgroups",
"",
"**INTERPRETATION AND REPORTING:**",
"☐ Score conversion and interpretation guidelines",
"☐ Reliability and precision estimates",
"☐ Limitation and assumption documentation",
"☐ Practical significance assessment"
)
for(check in implementation_checklist) {
cat(check, "\n")
}## **PRE-ANALYSIS CHECKS:**
## ☐ Adequate sample size (200+ respondents minimum)
## ☐ Sufficient category usage (5+ responses per category)
## ☐ Data quality (missing data patterns, outliers)
## ☐ Assumption verification (unidimensionality, local independence)
##
## **MODEL FITTING:**
## ☐ Software selection and parameterization
## ☐ Starting value specification
## ☐ Convergence monitoring
## ☐ Parameter estimate reasonableness
##
## **MODEL EVALUATION:**
## ☐ Overall model fit assessment
## ☐ Item-level fit diagnostics
## ☐ Person-level fit analysis
## ☐ Parameter stability across subgroups
##
## **INTERPRETATION AND REPORTING:**
## ☐ Score conversion and interpretation guidelines
## ☐ Reliability and precision estimates
## ☐ Limitation and assumption documentation
## ☐ Practical significance assessment##
## === COURSE COMPLETION ===cat("Bu dərsi tamamlayaraq, siz Graded Response Model və polytomous IRT-nin əsas prinsiplərini öyrəndiniz.\n")## Bu dərsi tamamlayaraq, siz Graded Response Model və polytomous IRT-nin əsas prinsiplərini öyrəndiniz.cat("İndi Likert scales və rating data-nı sophisticated IRT methods ilə analiz edə bilərsiniz.\n\n")## İndi Likert scales və rating data-nı sophisticated IRT methods ilə analiz edə bilərsiniz.## Next Steps for Polytomous IRT Mastery:## • Practice with ltm, mirt, and eRm packages## • Explore nominal response models for unordered categories## • Study multidimensional polytomous IRT models## • Learn about mixture polytomous IRT models## • Investigate computerized adaptive testing with polytomous items## • Apply polytomous IRT to real survey and assessment data## • Compare polytomous IRT with classical methods## Polytomous IRT opens new possibilities for sophisticated measurement!## Müvəffəqiyyətlər polytomous IRT expertise journey-nizdə!